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A two-dimensional lattice gas of two species, driven in opposite directions by an external force, 
undergoes a jamming transition if the filling fraction is sufficiently high. Using Monte Carlo simula- 
tions, we investigate the growth of these jams ( "clouds" ) , as the system approaches a non-equilibrium 
steady state from a disordered initial state. We monitor the dynamic structure factor S(k x , k y ;t) 
and find that the k x = component exhibits dynamic scaling, of the form S(0,k y ;t) = t^ S(k y t a ). 
Over a significant range of times, we observe excellent data collapse with a = 1/2 and /3 = 1. The 
effects of varying filling fraction and driving force are discussed. 

PACS numbers: 05.70.Ln, 68.43.Jk, 64.60.Cn 



I. INTRODUCTION 

The study of phase separation and coarsening in sys- 
tems undergoing continuous or first-order phase transi- 
tions has a long history in physics and materials science 
[1-4]. A model system, like an Ising lattice gas, or a real 
alloy, like a mixture of tin and lead, is prepared in a high- 
temperature state and then suddenly quenched below its 
co-existence curve. As the system phase separates, its 
properties are dominated by the morphology of growing 
single-phase domains. A particularly interesting feature 
of many phase-ordering systems is dynamic scaling: if 
space and time arc appropriately rcscalcd, growing do- 
mains at different times are found to be statistically self- 
similar, and the characteristic domain size, R(t), grows 
as a power of time, t a . More detailed information is con- 
tained in the equal-time two-point correlation function, 
or equivalently, the structure factor, which are general- 
ized homogeneous functions of space and time. For sys- 
tems evolving towards terminal equilibrium states, the 
domain growth exponent a and the scaling behavior of 
correlation functions are fundamentally well understood 
[3]. 

The situation is very different for many-body systems 
evolving towards non- equilibrium steady states (NESS). 
Maintained far from equilibrium by some external force, 
for example couplings to multiple energy or particle reser- 
voirs, these systems carry nonzero fluxes. As a result, 
their stationary distributions lie outside the Boltzmann- 
Gibbs framework and are known only for a few special 
cases. Yet, nonequilibrium systems occur frequently in 
nature, particularly in many biological contexts. Not 
surprisingly, they display much richer behaviors than sys- 
tems in thermal equilibrium [5, 6], including a variety of 
pattern forming instabilities and first-order phase transi- 
tions, controlled by the external drive rather than a tem- 
perature variable. However, rather little is known about 
coarsening phenomena in such systems. Given that the 
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underlying dynamics violates a very fundamental sym- 
metry of equilibrium systems, namely, detailed balance, 
it is not immediately obvious whether features such as 
dynamic scaling or power law growth will persist when 
systems evolve towards terminal states which fall into the 
NESS class. 

As a first step towards a better understanding of coars- 
ening in such systems, it is instructive to investigate a 
few simple models, in the hope that these will gener- 
ate insights from which a more general theory can be 
built. Looking for candidates which fall into the NESS 
class, which are well characterized in other sectors of their 
phase diagram and exhibit coarsening in some parameter 
regime, we are naturally led to driven diffusive systems 
[5, 7]. These systems involve one, or several, species of 
particles, diffusing on a lattice subject to a differential 
bias and short-range interactions. Both the prototype, 
first introduced [7] as a deceptively trivial modification 
of the Ising lattice gas, and its variants display many 
surprising and counterintuitive phenomena [5]. A par- 
ticularly interesting modification involves models with 
two particle species driven in opposite directions [8-10] 
where "jamming" transitions emerge from biased diffu- 
sion alone. 

Let us very briefly survey earlier studies of domain 
growth and dynamic scaling in driven diffusive systems. 
The prototype model, an Ising-likc lattice gas in which 
the particles are "charged" and driven by an external 
"electric" field E, sustains a nontrivial particle current 
on a fully periodic lattice. Still, the order-disorder tran- 
sition of the undriven system survives, separating a disor- 
dered phase from a low-temperature phase which phase- 
separates into high- and low-density strips, aligned with 
the drive. If the system is quenched from a typical high- 
temperature state into the phase-separated sector of the 
phase diagram, coarsening of single-phase domains oc- 
curs [11, 12]. Some interesting morphological discrepan- 
cies between simulation data and results from a contin- 
uum theory [11] were eventually resolved [13]. Turning 
to two-species models, the onset of jamming separates a 
homogeneous, high-current phase from a spatially inho- 
mogencous, low-current phase. As in the single species 
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(a) (b) (c) 

Figure 1: Snapshots of an 800 x 800 system, at (a) t = 1024, 
(b) 4096, and (c) 16384, in units of MCS. E = 10. Positive 
(negative) particles are black (white); holes are blue. 



case, the jams take the form of strips of high particle 
density, but these are now aligned transverse to the field 
direction. At the late stages of the approach to the steady- 
state, the system typically exhibits several strips which 
coarsen until only a single strip remains in the long-time 
limit. 

Earlier work on dynamic properties has mostly focused 
on these late stages. Since the strips are (on average) 
uniform in the transverse direction, they are quite well 
described by a set of mean-field equations, in one space 
dimension and time [9, 14]. If the excluded volume con- 
straint is enforced rigorously, so that the particles are 
not allowed to swap places, the strips coarsen logarith- 
mically slowly [14, 15]. Another group of studies inves- 
tigates systems where the microscopic dynamics is al- 
ready restricted to one [16, 17], or quasi-one [18, 19], 
dimension. For interesting behavior to occur, particle- 
particle ("charge") exchanges must be permitted, albeit 
with a small rate, compared to particle-hole exchanges. 
Provided the model parameters are chosen appropriately, 
compact particle clusters form easily, and coarsen until 
one a single large cluster remains. By virtue of the charge 
exchange process, power law growth dominates here. 

In this article, we present the first study of fully two- 
dimensional coarsening in a two-species model with a 
strict excluded volume constraint. Starting from an ini- 
tially disordered configuration, the system parameters 
(density, bias) are chosen so as to favor a jammed phase. 
Almost immediately, small "clouds" (Fig. 1) of locally 
jammed particles form. The larger clouds then grow, at 
the expense of the smaller ones, until a large cloud per- 
colates along the transverse direction, forming a strip. 
Eventually, several strips emerge and compete with one 
another, on much slower time scales. We focus on the 
multi-cloud regime, long before the late-stage strip coars- 
ening regime sets in. We monitor the equal-time struc- 
ture factor S(\s.,t), as a function of wave vector k and 
time t, averaged over initial conditions and system his- 
tories. A range of system sizes, densities, and E'-values 
are studied. Since the field selects a specific direction, 
the y-axis, it is not surprising that the structure factors 
are anisotropic, in k x and k y . More remarkably, we find 
that the system exhibits good dynamic scaling in k v and 
t, provided k x is fixed at k x = 0. Assuming the scaling 



form S(0,k y ,t) = t /3 S(k y /t a ), the scaling exponents are 
found to be a = 1/2 and (3=1. For nonzero values of 
k x , or in the full (k, t) domain, we have not been able to 
achieve good data collapse. 

This paper is organized as follows. We first present the 
model, a set of diagnostic observables, and some tech- 
nical details of the simulations. Next, we discuss our 
simulation results and evidence for dynamic scaling. We 
conclude with some comments and open questions. 



II. THE MODEL AND ITS OBSERVABLES 

Our model is defined on a two-dimensional square lat- 
tice of size L x x L y with fully periodic boundary condi- 
tions. Two species of particles, referred to as "positive" 
and "negative" , reside on the sites of the lattice, subject 
to an excluded volume constraint. Hence, a given config- 
uration of the system can be labelled by a set of occupa- 
tion variables, er(r), taking the values 0, +1, and —1 if 
the site r = (x,y) is empty or occupied by a positive or 
negative particle, respectively. The particles experience 
no interactions, apart from respecting an excluded vol- 
ume constraint. For simplicity, we restrict ourselves to 
systems which are neutral: X) r <T ( r ) = 0- For later refer- 
ence, we also define the particle (as opposed to charge) 
occupation n(r) via 

n(r) = \a(v)\ (1) 

so that the total particle density ( "mass" ) m is given by 

m = (L x L y y J2r n ( r )- 

In the absence of the driving force, the particles per- 
form simple diffusion, i.e., jump with equal probability to 
a randomly selected nearest-neighbor site, provided it is 
unoccupied. As a result, there is no net current (of cither 
mass or "charge") through the system, and the steady 
state is spatially uniform. In contrast, an "electric" field, 
applied in the positive y-direction, biases positive and 
negative particles in opposite directions. In our simula- 
tions, a bond is selected at random and the occupancies 
of the two associated sites are checked. If the bond carries 
a particle-hole pair, an exchange will always be made if 
this results in a positive (negative) particle moving in the 
transverse or positive (negative) y-dircction; otherwise, 
the exchange is attempted with rate exp(— E). Clearly, 
this dynamics is translation-invariant and invariant un- 
der charge-parity transformation (a — » —a, y — > —y). 
We use a random sequential dynamics, with one Monte 
Carlo step (MCS) corresponding to L x x L y update at- 
tempts. All runs start from a random initial condition. 

The system sizes studied ranged from 100 x 100 to 
3200x3200. The density varied from m = 0.3 torn = 0.7. 
We also considered different values for the probability for 
a particle to move backwards. Our reference system, for 
which the largest data set was collected, is an 800 x 800 
lattice, with m = 0.5 and E = 10. The latter gives 
a probability of 4.5 x 10 ~ 5 for backward jumps which 
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is zero for all practical purposes. Runs lasted at least 
8196 = 2 13 MCS, and data are typically averaged over 
1, 000 runs, unless stated otherwise. Time is measured in 
MCS. 

The final stationary state of the system is well under- 
stood. For sufficiently large particle density m and field 
E, the system displays a single strip of particles, trans- 
verse to the field direction [8]. In that fashion, trans- 
lational symmetry is spontaneously broken. The strip 
itself is charge-segregated, with positive (negative) par- 
ticles occupying sites with lower (higher) y-coordinatcs. 
The interior interface (separating positive from negative 
particles) is glassy, due to the absence of any charge ex- 
changes. In contrast, the exterior interfaces (separating 
particles from holes) is quite smooth, since its fluctua- 
tions are controlled by suppressed particle moves, i.e., 
by the parameter exp(— E). This parameter also con- 
trols the density of particles in the remainder of the 
system, reminiscent of a gas-liquid interface under grav- 
ity. Due to the periodic boundary conditions, however, 
a small current flows, even in the jammed phase, limited 
by e~xjp(—EL y ). A simple mean-field theory allows us to 
compute average density profiles, currents, and the phase 
diagram, in good agreement with the simulations [8-10]. 

In the following, we will always choose system param- 
eters such that the system evolves towards an inhomoge- 
neous, jammed final state. Starting from a random ini- 
tial disordered configuration, small jams of positive and 
negative particles form very rapidly, due to local density 
fluctuations. Some of these, typically the larger ones, will 
collect more particles and grow, while others shrink and 
dissolve, as illustrated in Fig. 1. We refer to this stage 
as the coarsening of clouds, or clusters. Eventually, first 
one and then several of the largest clouds will span the 
lattice in the transverse direction, and the evolution is 
no longer dominated by the coarsening of well-separated 
clouds. Now, multiple strips compete for particles until 
only a single one remains, and the system has reached its 
steady state. 

In this study, the characteristic shapes and separations 
of the clusters are of interest. Thanks to translational in- 
variance, a suitable observable is the equal-time structure 
factor, defined through the Fourier transform of the local 
occupation, 



S(k x , ky' : t) — 



L X Ly 



x=0 y=Q 



(2) 



where k x = 2t:1/L x , I = 0, 1, L x — 1, and k y = 2nj / L y , 
j = 0, 1, L y — 1. Here, n(x, y; t) denotes the local oc- 
cupation of site (x, y) at Monte Carlo time t. The aver- 
age (...) is taken over multiple runs, using configurations 
recorded at the same Monte Carlo time. All initial con- 
ditions are random. 

Let us first establish a few properties of this structure 



factor. The value at the origin is easily found: 
1 



S(0,0) 



L X Ly 



< \^n{x,y)\ 2 > = m 2 L x L y (3) 



Further, S is related to the two-point correlation func- 
tion, G{x,y) = (n(x,y;t)n(0,0]t)), via 



S(k x , k y ;t) = ^ G(x, U) ( 



i(k:cX+kyy) 



which also provides us with the sum rule 
^ ' S(k x , k y ;t) — mL x L y 



(4) 



(5) 



Finally, it is useful to evaluate S for a few special cases, 
including the initial and final configurations. Since the 
time argument is inessential here, it will be suppressed 
for now. If the system is filled randomly with particles, 
at density m, the structure factor is easily found to be 

S(k x ,k y ) = m(l-m)[l + 0(l/(L x L y ))] 

+ m 2 (L x L v )S k;cfi S ky:0 (6) 

Clearly, S(k x , k y ) is uniform for all non-zero k = {k x , ky). 

For comparison, we also evaluate the structure factor 
for a perfectly ordered single strip which reflects the sta- 
tionary state, modulo fluctuations: 



S(k x ,k y ) = y^4 x ,o 

Ly 



1 — cos (mkyLy) 



1 — cos k„ 



(7) 



These expressions provide a few benchmarks for the sim- 
ulation data presented below. 



III. SIMULATION RESULTS AND TESTS FOR 
DYNAMIC SCALING 

A. Unsealed structure factors. 

In this section, we first present Monte Carlo data for 
raw (unsealed) structure factors. We have collected data 
for a wide range of k. Roughly speaking, the k-value 
of the peak position reflects a characteristic spacing of 
the growing clusters, while the peak width carries in- 
formation about fluctuations. For illustration purposes, 
we show two projections here, namely S(0,k y ;t) and 
S{k x ,0; t) (Fig. 2). Plotted vs k y , the data for 5(0, k y ; t) 
show a distinct maximum which moves towards smaller 
values of k y for later times. In contrast, S(k x ,0;t) is 
monotonically decreasing in k x , for all t. A fuller pic- 
ture can be gleaned from contour plots of S(k x ,k y ;t) 
which indicate that, even for the earliest times consid- 
ered (t < 50 MCS), the maximum of S is found on the 
k x = axis. As time progresses, the peak position shifts 
from larger values of k y to smaller ones, and the peak 
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Figure 2: Unsealed structure factors S(0,k y ;t) (a) and 
S(k x ,Q;t) (b) for an 800 x 800 system, at E = 10. The dif- 
ferent curves correspond to different values of j (top) and I 
(bottom), specified in the legend. The notation is that of 
Eq. (2). Time is given in units of MCS. 



height increases. These findings suggest that, as the first 
clouds emerge from the fully disordered initial configura- 
tions, they quickly develop a characteristic length scale 
in the field direction, but remain disordered in the trans- 
verse direction. 



Scaled structure factors 



Snapshots of typical configurations at different times 
(Fig. 1) show clusters of particles ("clouds") which grow 
in both the parallel and the transverse directions. If a 
simple rescaling of system size renders configurations, 
recorded at different times, statistically similar, we can 
hope for dynamic scaling, as illustrated by Fig. 3. After 
an appropriate rescaling of Figs, lb and c, Fig. lc is plot- 
ted inside Fig. lb which is plotted inside Fig. la. One has 
to take a very careful look, if one wants to discern the 
internal boundaries (discontinuities) between the three 




Figure 3: Scaled configurations from Fig. 
details. 



See text for 



pictures. This illustrates - at a simple visual level - how 
closely they resemble one another, after rescaling. How- 
ever, our visual ability to detect scaling is easily deceived 
and provides, at best, the motivation for a more quanti- 
tative study. 

A quantitative test of dynamic scaling requires a care- 
ful analysis of the structure factors. Assuming that char- 
acteristic lengths in both directions increase as powers of 
time, but with possibly different exponents due to the 
anisotropy induced by the field, we first seek dynamic 
scaling in the form 



S(k x , ky]i) ~ t f{k x t 1 



! kyt 



(8) 



where the ~ indicates that we should expect this form to 
hold only for certain ranges of time and wavevector. The 
sum rule, Eq. (5), immediately leads to the exponent 
identity 



d\ + OL2 = (3 



(9) 



If dynamic scaling holds, one should be able to determine 
a set of scaling exponents in such a way that structure 
factor data for different times and wavevectors collapse 
onto a single curve if plotted according to Eq. (8). How- 
ever, we have not been able to achieve satisfactory data 
collapse for this general form. Once again this suggests 
that there are no characteristic transverse length scales, 
associated with this growth process. It also illustrates 
that merely visual tests of scaling, such as Fig. 3, must 
be treated with some caution. 

Turning to the remnant structures in S(0,k y ;t), the 
data in Fig. 2 show a sequence of curves of similar shapes, 
with the maximum shifting to smaller k y for later times. 
Even if the general form, Eq. (8), is not obeyed, we can 
explore the possibility of dynamic scaling in the reduced 
space k x = 0. In the remainder of this article, we focus 
on tests of 



S(0,k y ;t) ~t p f{0,kyt a ) 



(10) 



Fig. 4 shows the scaling plot for a half-filled system for 
times ranging from t = 2 9 = 512 to t = 2 13 = 8192. We 
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Figure 4: Scaling plot for S(0,k y ;t) for an 800 x 800 system, 
at m = 0.5 and E = 10. Five different times, ranging from 
t = 512 to t = 8192, in units of MCS, are shown, a = 1/2 
and (9 = 1. The solid line denotes a z -3 power law. 



find excellent data collapse with the scaling exponents 
a = 0.50 ± 0.02 and /3 = 1.00 ± 0.02. Much longer runs 
(with poorer statistics) show that the data continue to 
collapse well, until at least t ~ O(10 6 ). Our value for 
a, the exponent controlling the characteristic spacing of 
domains, stands in stark contrast to its counterpart for 
conserved coarsening in equilibrium systems. There, it 
takes the value 1/3, for a simple scalar density such as 
ours. 

The scaling function exhibits Gaussian behavior near 
the maximum, and falls off as z~ 3 , where z = k y t a is the 
scaling variable. This large z-behavior is highly reminis- 
cent of the Porod tail [20], well known in the theory of 
domain growth in equilibrium systems. There, it emerges 
from two essential features, namely, first, the presence of 
a single (isotropic) large length scale in the system, char- 
acterizing both the size and the separation of the coarsen- 
ing domains, and second, the existence of microscopically 
sharp domain walls. Here, the situation is more com- 
plex. While we do observe sharp domain walls between 
our clouds and the surrounding (nearly) empty regions, 
our model is manifestly not isotropic. What complicates 
the issue further is the absence of a characteristic spacing 
in the direction transverse to the field. Clearly, a more 
careful study is required before the large z-behavior of 
our model can be traced directly to a simple Porod law. 
As for the small z-behavior, we hesitate to offer any con- 
clusions. Certainly, it does not appear to follow the fc 4 
power law which would be expected for conserved coars- 
ening in equilibrium systems [21]. 

In the following, we probe the universality of the scal- 
ing exponents, as we change system parameters such as 
the particle density or the driving force. First, we con- 
sider the effect of system size. Since the scaling variable 



z defines a characteristic length scale £ tx t a , it is natural 
to expect a breakdown of scaling when £ becomes of the 
order of L, or when considering times t > L x l a . Indeed, 
scaling plots for a range of L, with 200 < L < 566 confirm 
this expectation very clearly. For example, in a 200 x 200 
system, the data for t = 10 12 = 4096 already deviate no- 
ticeably from the scaling curve, whereas for L = 566, 
such deviations are not observed until t = 2 15 = 32768. 

Next, we investigate the role of the overall particle den- 
sity. For coarsening in conserved equilibrium systems, it 
is well known that the scaling function depends on vol- 
ume fraction of the minority phase; however, the scal- 
ing exponents describing the structure factor remain un- 
changed [22]. Here, the situation is much more dramatic. 
Using the scaling exponents a = 1/2 and f3 = 1, the data 
collapse for densities close to half-filling (m = 0.40 and 
0.55) is still acceptable, but becomes progressively worse, 
for both larger (to = 0.70) and smaller (m — 0.30) densi- 
ties. Better data collapse can still be achieved, but at the 
price of modifying the scaling exponents. Fig. 5 shows 
the scaled data, with appropriately adjusted values of a 
and (3. It is natural to assume that these values reflect ef- 
fective, rather than true asymptotic exponents. A better 
understanding of the scaling function would be necessary 
to disentangle its TO-dcpcndence from the overall scaling 
exponents. 

We encounter a similar situation when considering the 
effect of the driving force, E. We find good data collapse, 
with a = 1/2 and (3 = 1, as long as the rate for a particle 
to move against its preferred direction, set by exp(— E), 
remains small. Once exp(— E) becomes comparable to 
0.2, deviations from scaling become noticeable. More 
work will be required to shed light on these preliminary 
observations. 



IV. CONCLUSIONS 



To summarize, we have explored the possibility of dy- 
namic scaling in a two-dimensional driven lattice gas, in- 
volving two species of particles. Positive and negative 
particles preferentially move in opposite directions and 
form small jams, due to an excluded volume constraint. 
Above a certain threshold density, these jams coarsen 
until a single strip of particles spans the system in the 
transverse direction. For an extended period of time, 
this coarsening process obeys dynamic scaling, provided 
we focus on characteristic length scales in the longitudi- 
nal direction. Monitoring a structure factor, S(0,k y ;t), 
we find very good data collapse provided t~@S is plotted 
vs k y t a . At and near half-filling (to = 0.5) and for large 
driving force, we find a = 0.50±0.02 and (3 = 1.00±0.02. 
For smaller E and densities further away from half-filling, 
we believe that the scaling function acquires a depen- 
dence on to and E. We note that we can still achieve rea- 
sonable data collapse with the simple form given above, 
but only at the price of adjusting the exponents a and [3. 
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Figure 5: Scaling plots for 5(0, k y ;t) for an 800 x 800 system, at E = 10. Four different densities are shown, and the scaling 
exponents are adjusted to give satisfactory data collapse: top left, m — 0.30, a = 0.46, ft = 0.98; top right, m — 0.40, a = 0.50, 
13 = 1.00; bottom left, m = 0.55, a = 0.50, /3 = 1.00; bottom right, m = 0.60, a = 0.44, (3 = 0.95. All errors are at most 5%. 
The solid lines denote z~ 3 power laws. 



We believe that these effective exponents mask possibly 
significant modifications to the scaling function. 

Naturally, a better analytic understanding of the ex- 
ponents and of the scaling function would be desirable. 
It will be interesting to see what future studies in both 
simulations and analytics would reveal. The observation 
that a is essentially 1 /2 points towards a diffusive mech- 
anism. Based on visual inspection alone, the clusters 
evolve by exchanging particles with one another. If this 
process is truly random - i.e., a cluster gains and loses 
particles with a fixed, constant rate, one should indeed 
expect to find a diffusive growth of characteristic length 
scales. Due to the drive, the particle exchange occurs pre- 
dominantly between clusters which are nearest neighbors 



in the transverse direction; hardly any interactions oc- 
cur between nearest neighbors in the transverse direction. 
This may explain the absence of any apparent structures 
in k x . Work is in progress to analyze a well-established 
mean-field theory for this model, in the hope of gaining 
a better understanding of exponents and scaling func- 
tion. If successful, it should also elucidate the deviations 
and similarities of our coarsening process with respect to 
those in equilibrium systems. 
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